Nonlinear dynamics of a shear banding interface 
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We study numerically the nonlinear dynamics of a shear banding interface in two dimensional 
planar shear flow, within the non-local Johnson Segalman model. Consistent with a recent linear 
stability analysis, we find that an initially flat interface is unstable with respect to small undulations 
for sufficiently small ratio of the interfacial width I to cell length Lx- The instability saturates 
in finite amplitude interfacial fluctuations. For decreasing l/L^ these undergo a non equilibrium 
transition from simple travelling interfacial waves with constant average wall stress, to periodically 
rippling waves with a periodic stress response. When multiple shear bands are present we find 
erratic interfacial dynamics and a stress response suggesting low dimensional chaos. 

PACS numbers: 47.50.+d, 47.20.-k, 36.20.-r. 



Complex fluids such as polymers, liquid crystals and 
surfactant solutions have mesoscopic structure that is 
readily perturbed by flow 0. For example, wormlike 
surfactant micelles with lengths of the order of microns 
can be induced to stretch, disentangle and entangle, and 
break or increase in length. Their mechanical response is 
therefore highly non-Newtonian, with shear flows induc- 
ing normal stresses {e.g. axx — '^yy), ^^nd with the shear 
stress Gxy being a nonlinear function of applied shear rate 
7. Recent work on such fluids has led to a fairly consis- 
tent picture of "shear banding" : coexistence in shear flow 
of viscously thicker (nascent) and thinner (flow-induced) 
bands of material flowing at different local shear rates, 
for an overall average imposed shear rate. 

This phenomenon can be described by constitutive 
models for which the shear stress is a non-monotonic 
function of shear rate for homogeneous flow. This leads 
to a separation into bands of differing shear rate that 
coexist at a common shear stress |2, 13 • Although 
most studies have assumed a flat interface between the 
bands and (with a few exceptions 0]) predicted a time- 
independent banded state, an accumulating body of data 
has demonstrated that the average shear str ess, and the 
banding interface, can fluctuate [1 13. II II ITa llll [l9l| . 
An important question is then whether these fluctuations 
resemble small amplitude capillary waves stabilised by 
surface tension, or whether they arise from an underly- 
ing instability, stabilised at large amplitude by nonlinear- 
ities. In this Letter we give strong evidence supporting 
the latter scenario, via the first theoretical study of the 
nonlinear dynamics of a shear banding interface. 

The model - The generalised Navier Stokes equation 
for a viscoelastic material in a Newtonian solvent of vis- 
cosity 77 and density p is: 

p(at-t-v.V)v = V.(S + 277D-FI), (1) 

where v(r) is the velocity field. The pressure P is deter- 



mined by incompressibility, V • v = 0. The viscoelastic 
stress S(r) evolves according to the non-local ("diffu- 
sive") Johnson Segalman (DJS) model [T^ IT^ 



(a* + V • V) S - a(D ■ S + S ■ D) 



(2) 



-(S ■ n - • S) = 2GD 



with plateau modulus G and relaxation time t. D and 
O are the symmetric and antisymmetric parts of the ve- 
locity gradient tensor, (Vv)^^ = daVp. For a — 1 and 
£ — Q this model reduces to the Oldroyd B model, which 
is motivated by considering an ensemble of beads paired 
by springs (simplified polymer chains) into dumbbells. 
Stress is generated as the flow deforms the dumbbells, 
and is relaxed on the timescale r for the springs to regain 
equilibrium length. To capture shear thinning the DJS 
model invokes a "slip parameter" a with \a\ < 1 to give 
non-affine dumbbell deformation The constitutive 

curve Txy ~ Y^xyij, a) +7?7 for homogeneous planar shear 
V = J/7X is then capable of non-monotonicity, allowing 
shear banding. The non-local diffusive term in Eqn.[2|ac- 
counts for spatial gradients across the banding interface 
on a length scale £ JJj. It arises naturally in models of 
liquid crystals, and diffusion of polymer molecules 1^. 
In the context of one dimensional (ID) calculations, it 
has been shown to give a unique (selected) stress T*y at 
which banding occurs |l3l |. as seen experimentally. 

We study 2D boundary-driven flow between parallel 
plates at y — 0, L with plate conditions dyT,ai3 = V a, /3 
for the viscoelastic stress, and no slip or permeation for 
the velocity. In the flow direction we consider a domain 
X — {0,Lx}, with periodic boundary conditions. We 
choose a = 0.3, rj = 0.05 throughout, and units in which 
G = l,r = 1 and L = 1. 

For an imposed average shear rate 7 = [vy{L) — 
Vy{0)]/L in the region of decrea sing stress, dTxy/dj < 0, 
homogeneous flow is unstable A ID (y) calcula- 
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tion then predicts separation into bands of shear rates 
7i = 0.66,72 — 7.09, at a selected shear stress T*y — 
0.506. Recent analysis showed this stationary ID 
banded state to be linearly unstable to 2D {x, y) per- 
turbations corresponding to undulations of the interface 
with wavevector q = qx^- The most unstable mode has 
QxL « 2n, and the instability involves feedback of the 
normal stress with velocity fluctuations across the inter- 
face. In this Letter, we study the fate of the interface 
in the nonlinear regime, and demonstrate it to be resta- 
bilised at the level of finite amplitude undulations. 
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FIG. 1: (a) Growth rate uj{qx) of perturbations about an 
initially ID banded state, (o, — ): predictions of linear sta- 
bility for f = 0.005, 0.01, 0.015. (□): early time dynam- 
ics of the full 2D numerics for e = 0.005, = 4, 6, 8; ^ = 
0.01, = 2, 4, 6, 8; ^ = 0.015, = 2, 4, 8. (b) Power spec- 
trum P(qx) in the travelling wave regime for I — 0.005, ~ 
1 (o); £ = 0.01, = 2 (□); e = 0.015, = 2, 4, 6 (A). 

Numerical scheme — We enforce incompressibility by 
eliminating the velocity in favour of a stream- function ip. 
At each time step we update the local parts of Eqn. [21 us- 
ing an explicit Euler algorithm within a finite difference 
scheme on a rectangular grid of {Nx,Ny) nodes, using 
third order upwinding for the convective term (v.V)S. 
We then take a Fourier transform x ya the flow 

direction and update the nonlocal part of Eqn. |5] using 
the semi-implicit Crank-Nicolson algorithm. We finally 
update Eqn. at zero Reynolds number, p = 0. We 
attain convergence to one percent with respect to in- 
creasing spatial and temporal resolution. This ensures 
that the results shown are converged to the eye, apart 
from the stress signal Txy{t) in Fig. 01 this shows slight 
quantitative, but not qualitative, changes. Runs at finer 
resolution are prohibitively time consuming. We study 
two initial conditions: (ICl) the stationary banded state, 
predicted by a ID calculation, of two bands separated by 
flat interface; and (IC2) a homogeneous unstressed fluid, 
corresponding to shear startup from rest. 

Linear regime — First we consider the early time evolu- 
tion of an initially flat interface, ICl. For a wide enough 
interface i and small enough system size this ID solu- 
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FIG. 2: "Phase diagram" for 7 — 2.0, showing the different 
non-linear regimes for ICl (accurate to spacing between num- 
bers) and the numbers of linearly unstable modes (Ref. .16]) ■ 



tion is stable. For smaller £/Lx, however, interfacial un- 
dulations are predicted by the linear analysis of Refs. 
to become unstable. Our numerics successfully repro- 
duce this instability during the initial evolution away 
from ICl: the eigenvectors and gro wth rates Lo{qx) match 
the analytical results of Refs. FigCJ Beyond this 
early-time regime, the interface is restabilised by nonlin- 
ear effects, and the instability saturates in undulations of 
finite large amplitude. Depending on the distance from 
the onset of linear instability. Fig. [3 the ultimate attrac- 
tor corresponds either to a steady travelling interfacial 
wave, or to periodically rippling waves. 

Travelling wave — For values of £/Lx marginally in- 
side the unstable regime (Fig. 121 the ultimate attrac- 
tor comprises a travelling wave A A(y, x — ct) for 
all order parameters A — {Sxx, ^xy, "^yyiip}- Since this 
satisfies periodic boundary conditions, the wall-averaged 
shear stress is constant in time, with a value Txy,ss — 
0.51 — 0.54 that depends on £ and Lx and is slightly 
higher than the selected stress T*, = 0.506 of the ID 
calculation. In Fig.ylthe wave-speed c = 1.15, compared 
with a horizontal velocity Vx that varies between 1.439 
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FIG. 3: Greyscale of order parameters for travelling wave in 
the [x, y) plane for £ = 0.015, Lx ~ 6, and upper wall velocity 
V = 7L = 2 to the right. 
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FIG. 4: Rippling wave at £ = 0.005, = 4, 7 = 2. Top: 
Greyscale of T^^cxix, y). Upper wall moves to the right. White 
line in d): interface height defined in text. Bottom: aver- 
age wall stress. Inset: {Nx,Ny,dt) = (200,800,0.0003), solid 
line; (400,800,0.0001), dotted; (400,800,0.0003), dashed; 
(500,800,0.0002), dot-dashed. 



and 1.484 along the interface. To analyse the structure 
of this state we define the following power spectrum, 



4 i-i 



dy\A^{q^,y)\ 



(3) 



As seen in Figs. ^ and |21 while the linear stability anal- 
ysis predicted only a few unstable modes, the ultimate 
nonlinear state has active modes at higher qx, Fig. 
Interestingly, the dominant mode is still close to that of 




FIG. 5: Rippling wave regime. Greyscale of interface height 
(defined in text) in {x,t) plane. Parameters as in Fig. |1] 
Time t — 150 . . . 160 upwards. Upper: raw data. Lower: 
transformed as x ^ x — ct, c extracted by eye. 



the linear analysis, despite being in a nonlinear regime 
with a finite interfacial displacement Sh. Indeed, for all 
cases in Fig. QJj, we find empirically that the dominant 
mode has the longest wavelength that is both consistent 
with periodic boundary conditions and linearly unstable. 

Rippling wave — For smaller £/Lx, deeper inside the 
unstable regime (Fig. [SJ, we see a new regime in which 
the travelling wave now periodically "ripples" . The wall- 
averaged stress Txy is periodic in time (Fig. ^ again for 
ICl), with variations of the order of one percent, and an 
average value larger than the ID selected stress T*y. The 
interface height h{x) is shown as a white line in Fig.QJi 
and as a greyscale over a time window in Fig. \S\ Because 
the structure is quite complicated near a rippling wave 
(Figs.0]:d), the location of h{x) depends on its precise 
definition. We take h(x) as that value of y at which "E^x 
lies half way between its values at the two walls. (An 
alternative might be J dy y \dyT,xxiy)\-) This provides a 
fairly reliable measure, subject to a small kink in the 
rippling region (Fig. ^) . 

Robustness to initial conditions — So far, we have con- 
sidered only ICl, two bands separated by a single flat in- 
terface. To test robustness to ICs, we now study startup 
from rest (IC2). Here the system develops either (i) 
two bands that show the same dynamics as with ICl, 
or (ii) multiple bands that have erratic dynamics sug- 
gestive of low dimensional chaos. Fig. El This is the 
counterpart of startup from rest in ID planar shear flow, 
which typically yields random configurations of multiple 
bands, because the uniform shear stress Txy{y) = T*y 
allows interfaces to reside at any y value. (To achieve 
ICl, the system was prebiased to form two bands.) In 
contrast, Couette flow between concentric cylinders has 
a stress T^y ~ (r=radius), which allows only a 

single stationary interface at r* 1/ y/T*y, as demon- 
strated in ID in Ref. [Hill. By analo gy, we anticipate 
that curvature should eliminate the multiple band case 
(ii) in 2D. Although a true Couette calculation has not 
been performed, we implement a "poor man's" version by 
adding a biasing stress-gradient in the y direction while 
the bands are forming. For strong enough bias, we then 
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FIG. 6: Multiple bands after startup from rest for i — 0.005, 
Lx = 2, -y = 2. Top: greyscale of 'Exxix,!/) at time t — 120, 
reconstructed from the first 15 Fourier modes. Bottom: Shear 
stress Txy{t). 

indeed find case (i), just two bands. An open question is 
whether, even with ICl, a transition exists in the numer- 
ically inaccessible regime of small i/L^ (far bottom right 
in Fig. 12) to chaotic dynamics of a single interface. 

To summarise, we have studied numerically the nonlin- 
ear dynamics of a shear banding interface in 2D planar 
shear flow, within the non-local DJS model, for a sys- 
tem of length with periodic boundary conditions in 



X. Upon decreasing the ratio of the interfacial width £ to 
the system length the stable two band state separated 
by a single interface undergoes successive transitions to 
travelling waves and rippling travelling waves. Multiple 
shear bands can also be found, depending on initial con- 
ditions, which show irregular interfacial dynamics and a 
corresponding stress signal suggestive of low dimensional 
chaos. Depending on the regime, this dynamics is quali- 
tatively similar to several recent experiments in wormlike 
micellar solutions that studied the interfacial dynamics 
and the associated stress response 0, 0, 0, . 

In earlier work. Yuan et al. [2^ used a 2D finite ele- 
ment algorithm to evolve a related JS model in planar 
shear. Contrary to our results, they reported a stable in- 
terface. However, they considered rather small values of 
Lx, perhaps in the stable regime, and the finite element 
technique may have introduced stabilising numerical dif- 
fusion. Related recent work by Onuki j^J on a two- 
fluid version of the DJS model |22| at imposed stress also 
showed a variety of interesting time-dependent behaviour 
due to the interplay between stress and concentration de- 
grees of freedom; it would be interesting to study the 
phenomena found here in that model. Indeed, an open 
question is the extent to which the instability seen in this 
DJS model is ubiquitous among other models of shear 
banding, such as those for liquid crystals. Future work 
should also extend to larger systems and smaller inter- 
faces, likely to lead to more complex dynamical states. 
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